Effect of bond length fluctuations on crystal nucleation of hard bead chains 
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We study the effect of bond length fluctuations on the nucleation rate and crystal morphology 
of linear and cychc chains of flexibly connected hard spheres using extensive molecular dynamics 
simulations. For bond length fluctuations as small as a tenth of the bead diameter, the relaxation 
and crystallization resemble that of disconnected spheres. For shorter bond lengths and chains with 
< 10 beads the nucleation rates depend sensitively on bond length, number of beads per chain, and 
chain topology, while for longer chains the nucleation rates are rather independent of chain length. 
Surprisingly, we find a nice data collapse of the nucleation rate as a function of number of bonds per 
sphere in the system, independent of bead chain composition, chain length and topology. Hence, the 
crystal nucleation rate of bead chains can be enhanced by adding monomers to the system. We also 
find that the resulting crystal morphologies of bead chains with bond fluctuations resemble closely 
those for hard spheres including structures with five-fold symmetry. 

PACS numbers: 82.70.Dd, 82.60.Nh, 64.70.D 



Although crystal nucleation from a supersaturated 
fluid is one of the most fundamental processes during 
solidification, the mechanism is still far from being well 
understood. Even in a relatively simple model system 
of pure hard spheres, the nucleation rates obtained from 
Monte Carlo (MC) simulations using the umbrella sam- 
pling technique difi'er by more than 6 orders of magnitude 
from those measured in experiments 111]. This discrep- 
ancy in the nucleation rates has led to intense ongoing de- 
bates in the past decade on the reliability of various tech- 
niques as employed in simulations and experiments to 
obtain the nucleation rates 0, Q . Recently, it was shown 
that the theoretical prediction of the nucleation rates for 
hard spheres is consistent for three widely used rare-event 
techniques, despite the fact that the methods treat the 
dynamics very differently Moreover, the structure of 
the resulting crystal nuclei as obtained from the different 
simulation techniques all agreed and showed that the nu- 
clei consist of approximately 80% face-centered-cubic-like 
particles. The predominance of face-centered-cubic-like 
particles in the critical nuclei is unexpected, as the free 
energy difference between the bulk face-centered-cubic 
(fee) and hexagonal-close-packed (hep) phases is about 
0.001 kgT per particle, and one would thus expect to find 
a random-hexagonal-close-packed (rhcp) crystal phase 
More surprisingly, simulation studies showed that 
the subsequent growth of these critical nuclei resulted 
in a range of crystal morphologies with a predominance 
of multiple twinned structures exhibiting in some cases 
structures with a five- fold symmetry [1, y]. Such struc- 
tures are intriguing as the fivefold symmetry is incom- 
patible with space-filling periodic crystals. Moreover, the 
formation mechanism of these fivefold structures is still 
unknown. Bagley speculated that the fivefold structures 
are due to the growth of fivefold local structures (a dec- 
ahedral or pentagonal dipyramid cluster of spheres) Jl\ 
that are already present in the supersaturated fluid phase 



[8|. Another mechanism that has been proposed is that 
these multiple twinned structures with a fivefold symme- 
try originates from nucleated fee domains that are bound 
together by stacking faults For instance, five tetra- 
hedral fee domains can form a cyclic multiple twinned 
structure with a pentagonal pyramid shape. A recent 
event-driven Molecular Dynamics simulation study on 
hard spheres showed, however, no correlation between 
the fivefold local clusters that are already present in the 
supersaturated fluid and the multiple twinned structures 
in the flnal crystal phase @. Hence, it was concluded 
that crystalline phases with multiple stacking directions 
may possess fivefold structures, whereas crystals with a 
unique stacking direction do not show any five-fold sym- 
metry patterns. These authors also succeeded to study 
the crystallization of tangent hard-sphere chains by in- 
troducing several clever (unphysical) chain-connectivity 
altering MC moves Surprisingly, they observed that 
tangent hard-sphere chains (no bond fiuctuations) never 
formed crystalline structures with a fivefold symmetry as 
the chain connectivity prohibits the formation of twinned 
structures and forces the crystals to grow in a single 
stacking direction d, [l3| ■ 

In this Letter, we study the effect of bond length fluc- 
tuations on the nucleation rate and crystal morphology 
of flexibly connected hard spheres. These hard-sphere 
chains with bond fluctuations can serve as a simple model 
for granular ball chains 11], colloidal bead chains [l3| 
or (short) polymeric systems. A better understanding 
of the behavior of these bead chains may shed light on 
the glass transition and crystallization of polymers. In 
fact, it has been shown recently that random packings of 
granular ball chains show striking similarities with the 
glass transition in polymers [ll]. Moreover, the scal- 
ing dependence of a number of physical characteristics 
on packing densities at concentrations which still remain 
unreachable in simulations of atomistic or coarse-grained 
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polymer models has been studied by employing this poly- 
meric hard-sphere chain model [isl . We also mention 
that recently, a colloidal model system of bead chains 
has been realized consisting of colloidal spheres that are 



where (t) is the average waiting time before a nucleation 
event occurs in a system of volume V. In order to iden- 
tify the crystalline clusters in the fluid phase, we employ 
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bound together with "flexible linkers" 

We consider a system of M bead chains consisting of 
N identical hard spheres with diameter cr in a volume V. 
The hard-sphere beads are connected by flexible bonds 
with a bond energy Ubond{i'ij) given by 



Ubondjrij) 

knT 



a < Tij < a + 6 
oo otherwise 



(1) 



where is the center-to-center distance between two 
connected spherical beads i and j, S is the maximum 
bond length, the Boltzmann constant, and T the tem- 
perature. The maximum bond length S varies from to 
O.lcr in our simulations, such that (5 = corresponds 
to a freely jointed chain of tangent hard spheres. We 
employ event driven molecular dynamics (EDMD) sim- 
ulations [3, [3 instead of the extensive MC method [§| 
to mimic the dynamics of granular and colloidal bead 
chains with flexible bonds. Since the pair potentials be- 
tween all spheres are discontinuous, the pair interactions 
only change when the beads collide or when the maxi- 
mum bond length is reached. Hence, the particles move 
in straight lines (ballistically) until they encounter an- 
other particle or reach the maximum bond length. The 
particles then perform an elastic collision. These colli- 
sions are identified and handled in order of occurrence 
using an EDMD simulation JJJ. 

We study the effect of bond length fluctuations on the 
nucleation of bead chains at a packing fraction rj — 0.55, 
which is at the upper boundary of the nucleation regime 
of hard spheres. At 77 = 0.55, the nucleation time ex- 
ceeds the relaxation time, which enables us to equili- 
brate the metastable fluid and to measure the nucleation 
rate. In order to obtain the initial configuration for our 
EDMD simulations, we use the Lubachevsky-Stillinger 
algorithm [isj to grow the beads of the chains in the 
simulation box to the packing fraction of interest with a 
very fast speed, i.e., da{t)/dt = O.Ola/r where a{t) is 
the size of spheres at time t with a the target sphere size 
and T = a^Jm/kBT the MD time scale. We monitor the 
relaxation of the initial configurations by comparing the 
bond angle distributions with those of previous work [l9| 
and by calculating the long-time self diffusion coefficient 
Dl of the beads. We find that Dl increases linearly with 
1/A^ as predicted by the Rouse model [2^, and is rather 
independent of the bond length 5, at least for the range 
of values that we studied 21]. In order to exclude the 
effect of dynamics on the nucleation rates 0, H^l , we use 
the long-time diffusion time of the beads, tl = 1/6Dl, as 
the unit of time in the nucleation rate. Starting from the 
metastable fluid phase, an EDMD simulation is employed 
to evolve the system until a spontaneous nucleation event 
occurs. The nucleation rate is then given by / = l/{t)V, 
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the local bond-order parameter analysis 

We calculate the crystal nucleation rates in systems 
consisting of ~ 10** beads for various linear and ring-like 
hard-sphere polymers for chain length 1 < < 20 and 
various bond length d. We mention here that we did 
not exclude inter-ring catenations in the configurations 
of ring polymers, which may occur for rings with N > 5. 
Figure [1] shows that for 6 as small as O.lcr, the nucleation 
rate resembles that of loose hard spheres. For 6 < 0.1a, 
we find that the nucleation rate decreases by several or- 
ders of magnitude with decreasing 6 and increasing chain 
length TV, but remains constant for sufficiently long poly- 
mers, i.e., N > 10. To check this surprising result, we 
also determined the nucleation rate of a single polymer 
of length N = 10'' and 5 = O.OScr, where all the beads are 
doubly connected except the two end beads. We observe 
that the system remains in the fiuid phase for ~ 7000t, 
before a critical nucleus of ~ 100 beads forms in the mid- 
dle of the chain, which subsequently grows further until 
the whole system is crystalline. In Fig. [21 we show the 
size of the largest crystalline cluster as a function of sim- 
ulation time from a typical MD trajectory. Additionally, 
we find that the nucleation rate does not decrease signif- 
icantly for A^ = lO"', which is highly unexpected as the 
beads can only move collectively. 

Additionally, we did not observe spontaneous nucle- 
ation within our simulation times for linear bead chains 
with 6 = 0.04cr and N > 5, which is to be expected 
as the bead chains become too frustrated on an fee lat- 
tice at ?7 = 0.55. At this density, the lattice spacing is 
much larger than the sum of the maximum bond length 
and bead diameter and thus the nucleation is strongly 
inhibited by geometrical frustration. However, for higher 
packing fraction (lattice spacing similar to 5 -I- cr), the 
nucleation barrier vanishes at such a high supersatura- 
tion and crystallization sets in immediately: a nucleation 
rate is here undefined. For lower 77, the nucleation times 
for bead chains exceed rapidly our simulation times, and 
again the rate cannot be determined. Finally, we also de- 
termine the nucleation rates for cyclic bead-chains (ring 
polymers) with bond length S = 0.05cr. Figure [1] shows 
that the nucleation rate of ring polymers does not de- 
crease monotonically with A'^ and is always lower than 
for linear polymers with the same length. However, the 
difference in nucleation rate is small for A^ > 10. 

As mentioned above, our results show that the nucle- 
ation rates of linear polymers decreases with chain length 
A^. One may argue that the nucleation rate is largely de- 
termined by the chain connectivity or the average num- 
ber of bonds per sphere in the system. In Fig. 3, we plot 
the nucleation rate for linear bead chains with maximum 
bond length 5 = 0.04cr and chain length N = 1,2,3,4, 
and 5, which correspond to an average number of bonds 
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FIG. 1. (Color online) Nucleation rates la^TL of linear hard- 
sphere polymers with maximum bond length S — 0.1a (green 
triangles), 0.05(7 (solid circles), 0.04a (solid squares), and of 
ring-like hard-sphere polymers with maximum bond length 
5 = 0.05a and chain length (open squares). For comparison, 
we also plot the nucleation rate for hard dumbbells with a 
maximum bond length S — 0.02a (blue triangles) and 5 = 
(diamonds). 
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FIG. 2. (Color online) Size of the largest crystalline clus- 
ter rimax as a function of simulation time for a single linear 
hard-sphere polymer with chain length A*' = 10*. The insets 
are snapshots at f = 7060t, SOOOr and lOOOOr, respectively, 
where only the solid-like beads are shown, and the two dark 
yellow spheres are the two ends of the chain. 



per bead of nt = {N - l)/N = 0, 1/2, 2/3, 3/4 and 4/5. 
In order to investigate the effect of average number of 
bonds per bead in the system on the nucleation rate, 
we also perform simulations for binary mixtures of linear 
and ring polymers with different chain lengths and bond 
length 5 = 0.04cr. We consider mixtures of linear chains 
with length TV = 2 and 5, TV = 1 and 6, and iV = 1 
and 10, and a mixture of monomers and ring polymers 
with = 10 and 20. The composition of the mixture is 
chosen such that the value of nt matches with one of the 
values for the pure systems. We compare the nucleation 
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FIG. 3. (Color online) Nucleation rates la^TL of linear hard- 
sphere polymers with chain length A'^ = 1,2,3, and 5 and 
of binary mixtures of linear hard-sphere polymers with chain 
length N = 2 and 5, A'^ = 1 and 6, and A'' = 1 and 10, as a 
function of the average number of bonds per bead rib- The 
composition of the mixture was chosen such that the value for 
nj, matches with one of the values for the pure systems. The 
maximum bond length equals S = 0.04a for all bead chains. 
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FIG. 4. (Color online) Typical configurations of the crystal 
structures for linear hard-sphere chains with chain length A'' = 
20 (a) and for ring-like polymers with A' = 3(b) and A' = 5(c). 
Only crystalline spheres are shown here. The blue and red 
spheres are hcp-like and fcc-like particles, respectively. 



rates for the pure and binary systems in Fig. |3l and ob- 
serve a nice data collapse, suggesting that the nucleation 
rate is determined by the number of bonds per sphere 
in the system, and the frustration imposed by the chain 
connectivity in these systems. We wish to remark here 
that no nucleation was observed in pure systems of lin- 
ear or ring polymers with 6 = 0.04cr and A'^ > 10 within 
our simulation times. Thus the addition of monomers 
enhances significantly the nucleation of polymers. 

Additionally, we investigate the structure of the result- 
ing critical nuclei and the full-grown crystals by calculat- 
ing the averaged local bond order parameters and W4 
for each sphere i that has Nb{i) > 10 neighbours. This 
analysis allows us to check whether a bead is fcc-like or 
hcp-like 



22l . [24|. We find that the critical nuclei con- 



tain more fcc-like than hcp-like particles, which is very 
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similar to the critical nuclei observed in hard-sphere Q 

In addition, we find 



22| 



and hard-dumbbell nucleation 
that the full-grown crystals display a range of crystal 
morphologies including random-stacked hexagonal close- 
packed crystals, and structures with a five-fold symme- 
try pattern for all polymer systems that we considered, 
even for ring-like polymers with chain length as small as 
iV = 3, which were not observed in the crystal of freely 
jointed chains of tangent hard spheres [6(. Exemplar- 
ily. Fig. |4] shows typical configurations of these five-fold 
symmetry patterns formed by linear hard-sphere chains 
of length = 20 and S = O.OScr, and cyclic bead chains 
of iV = 3 and 5 and S = 0.05(t. We note that the crystal 
structures resemble closely those observed in MD simula- 
tions of hard spheres [1, lij . As the crystal morphology is 
mainly determined by the crystallization kinetics rather 
than the bulk and surface contributions to the free en- 
ergy of the nucleus, it is tempting to speculate that the 
crystallization dynamics of hard bead chains with bond 
length fluctuations is similar to that of hard spheres, 
which is different from tangent hard-sphere chains (no 
bond length fluctuations) , which crystallize into stacking 
faulted layered crystals with unique stacking directions 
and no five-fold symmetry patterns 

Furthermore, we determine the bond angle distribu- 
tion function p(cos 0) in order to quantify the distribu- 
tion of bond angles between two neighboring polymer 
bonds in the supersaturated fluid and crystal nuclei. Fig- 
ure [5] shows p(cos 0) for systems consisting of linear and 
cyclic bead chains with N = 10. For the fluid phase, 
p{cos9) displays two peaks at cosd = —0.5 and 0.5, i.e., 
6 — 120 and 60 degrees, which corresponds with the 
most frequent three particle structures observed in ran- 
dom packings of spheres However, ^(cos^) of the 
crystal structures exhibits four pronounced peaks located 
around cos9 — —0.5,0,0.5,1.0 and a smaller peak at 
cose = V3/2 , which corresponds with 9 = 120,90,60,0 
and 30 degrees, respectively [19)]. In order to compare 
p{cos9) with that for a self-avoiding random walk on 
an fee crystal lattice, we integrate the four pronounced 
peaks of p(cos 0) between the two neighboring local min- 
ima. The results are shown in Fig. [5] together with the 
bond angle distribution function for a self-avoiding ran- 
dom walk on an fee crystal lattice. We find that p(cos 9) 
for crystals of linear and cyclic polymers agree well with 
that of a self-avoiding random walk Hence, p(cos 9) 
seems not to be affected by the chain connectivity of the 
polymer chains. Free energy calculations show indeed 
that the stable solid of freely jointed hard-sphere chains 
is an aperiodic crystal phase, where the spheres are posi- 
tioned on an fee lattice with the bonds randomly oriented 

In conclusion, we studied the effect of bond length fiuc- 
tuations on the nucleation rate and crystal morphology of 
linear and cyclic hard-sphere chains by varying the max- 
imum bond length parameter S. For sufficiently large S, 



3p- 
^2.5- 

0.5- 



0- 

,0.4- 



-0.3- 



Crystal of Ring polymer 

Fluid of ring polymer 

Crystal of linear polymer 

Fluid of linear polymer 



8 0.2 - 



0.1- 




e 



-0.5 





cos(0) 



I 

FCC 

Linear Polytner 
Ring Polymer 

\ , 

0.5 



FIG. 5. (Color online) Bond angle distribution function 
p(cos^) for a supersaturated fluid and crystal of linear and 
cyclic hard-sphere polymers with chain length = 10 and 
S — 0.05(j(top). Integration of the peaks of p{cos9) between 
the two neighboring local minima (bottom). For compari- 
son, we also plot the bond angle distribution for self-avoiding 
random walks consisting of 10 steps on an fee crystal lattice 
(bottom). The inset shows the definition of the angle 9. 



i.e. ^ 0.1(7, the relaxation and crystallization resemble 
that of disconnected spheres. Hence, we find nucleation 
rates and crystal morphologies similar to that observed 
for hard spheres. For bond length values < 6 < O.lcr, 
our results show on the one hand that the nucleation 
rates decreases significantly with S, but on the other 
hand, the crystal morphology of bead chains with floppy 
bonds resemble closely that of loose hard spheres 5, 6]. 
Apparently, the bond constraints slows down the relax- 
ation dynamics as the crystallization is frustrated by 
the chain connectivity, but docs not alter the crystal- 
lization kinetics. We thus flnd that the resulting crys- 
tal nuclei of bead chains with bond fluctuations pos- 
sess a range of crystal morphologies including structures 
with a flve-fold symmetry similar to hard-sphere crystals, 
while tangent hard-sphere chains (no bond length fluctu- 
ations) crystallize into stacking faulted layered crystals 
with unique stac king directions and no flve-fold symme- 
try patterns [13, Il9| Hence, the crystallization kinetic 
pathways of bead chains depends strongly on bond length 
fluctuations. Additionally, we flnd a nice data collapse of 
the nucleation rate as a function of number of bonds per 
sphere in the system, independent of bead chain compo- 
sition, chain length and topology. The crystal nucleation 
rate can thus be enhanced by the addition of monomers. 
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